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ABSTRACT 

The  Wiener  integral  formulation  combined  with  Monte  Carlo 
sampling  has  "been  used  to  compute  the  two-particle  Slater  sum  for 
He.  for  temperatures  ranging  from  273  K  down  to  2  K,  the  lower 
practical  limit  for  this  computational  method.   This  is  equivalent 
to  a  calculation  of  the  density  independent  part  of  the  pair 
distribution  function.   A  Lennard-Jones  6-12  potential  has  been 
used  to  describe  the  interaction.   Contributions  from  exchange 
were  found  negligible  at  5  K  and  above.   Comparisons  with  the 
Wigner-Kirkwood  expansion  are  made.   The  second  virial  coefficients 
derived  from  these  results  are  within  two  per  cent  or  three  per 
cent  of  the  results  obtained  from  the  usual  phase  shift  calculation. 
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1.   INTRODUCTION 


For  a  system  of  N  identical  particles  of  mass  m  enclosed  in  a 
volume  Q.,   with  Hamiltonian  H  ,  the  Slater  sum  is 


_OTT 

W..  =  N!  \3N  I  ¥*  (1,  2,    ...,  N)e   N  ¥.  (l,  2,  ...,  N),        (l.l) 

J^  ]_      <-0      «*Srf  *Nrf  X 


where  Y  (l,  ....  N)  is  the  wave  function  of  the  system  in  the  state  i;   1  is 
the  position  coordinate  of  particle  1,  2  is  the  position  coordinate  of 
particle  2,    etc.;  \   is  the  thermal  wave  length 

1 
,=  (^&)2;  (1.2) 


and 


e-fe-  (1'3) 


The  wave  functions  are  normalized  to  1  in  the  volume  Q,, 


y*   (1,    2,    ...,    N)  Y.  (l,  2,  . ...  N)d3ld32  ...  d3N  =  1,       (1.4) 


and  the  summation  in  Eq.  (l.l)  extends  over  all  states  appropriate  to  the 
statistics  of  the  system.   A  superscript  is  used  on  W  to  explicitly  denote 

"R  ~W 

Bose-Einstein  (W  )  or  Fermi-Dirac  (W  )  statistics. 


■p 

We  present  here  the  results  of  computing  W  for  ten  temperatures 
extending  from  273  K  down  to  2  K.   The  potential  describing  the  interaction 
is  the  Lennard-Jones  6-12  potential 


V  =  ta((2)12-  (^)6J,  ■   (1.5) 

2 
where  OC   and  a  are  the  deBoer,  Mich  els   values  appropriate  to  He.  : 


a  =   lU.OU  X  10~   erg, 


-ft 
a  =  2.56  X  10~  cm. 


(1.6) 


A  central  feature  of  this  calculation  is  the  use  of  the  Wiener  integral 

3 
formulation  of  the  Slater  sum,   described  in  the  following  section.   The 

Wiener  integrals  have  been  evaluated  by  a  Monte  Carlo  sampling  scheme  on  the 

ILLIAC  II  computer. 

The  results  exhibited  here  are  appropriate  to  a  description  of  the 

pair  distribution  function  n_(l,  2)  at  very  low  densities.   This  function  is 

given  by 

n2(i'  s)  = : — rw~fvn{~'  -'  ■■■'  ~)d3^'  d3~'  ••"  d3~'    {1-T) 

where  Q^  is  the  partition  function 

NA   a 


With  the  normalization  used  here 


,,   „n   N(N-l)    2 


*>#    *v 


(1.9) 


for 


l-2-oo, 


'  *w     *v 


(1.10) 


where  v  is  the  density.   At  suff iciently  low  densities  an  expansion  of  the 

k 
pair  distribution  function  in  powers  of  the  activity ,    z,    can  be  made: 


n(l,  2)  =  \^     Z   ^  (1,  2)z^+\ 

i=l 


(1.11) 


where  we  use  the  following  definition  of  z 


Q 


z  = 


N-l 


Q 


(1.12) 


'N 


The  functions  b„  are  modified  cluster  integrals;  for  £=1   and 


i=2  they    are   given  by 


b    (1,    2)    =  W    (1,    2), 

l_    *\*       «>^  cL    ^*       **+ 


(1.13) 


b2(i'  -}  =  ^3  {  /[  W3(^  ~>  $   -  W2(i>  S)W1(^)  ]  d^  }•     (l'lM 
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At  very  low  densities  the  activity  is  approximated  by 

z  ~  v  *.3  (1.15) 

and,  using  just  the  first  term  in  Eq.  (l.ll), 

I  n2(i'  &   *  *   W2(i'  &'  (l'l6) 

Because  of  spherical  symmetry  in  the  interaction,  Eq.  (1.5), 
W_(l,  2)  depends  only  on 

s  =  ll  -  2\,  (1.17) 


1  *\*     *v  ' 


and  so  we  may  write 

Wg(l,  2)  =  W2(S).  (1.18) 

It  is  tacitly  assumed  here  that  0,   is  so  large  that  boundary  effects  can  be 
ignored.   The  radial  distribution  function,  g(s),  in  the  approximation 
represented  in  Eq.  (l.l6),  is  given  by 

g(s)  sj  w2(s).  (1.19) 

Our  results  can  also  be  used  to  compute  the  second  virial 
coefficient,  given  by 

00 

B  =  -2ttN0J   (W2(S)  -  l)S2dS,  (1.20) 
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where  Nn  is  Avogadro's  number.   We  have  made  this  calculation  and  found 
good  agreement  with  other  calculations  of  B. 

A  secondary  reason  for  presenting  this  work  is  to  illustrate  the 
use  of  Wiener  integrals  as  a  computational  tool.   Although  the  Wiener  integral 
formulation  has  been  known  for  some  time  it  has  found  relatively  little  use 
as  a  computational  tool.   The  reasons  are  probably  two-fold;  although  it  has 
been  known,  it  has  not  been  well  known,  and  the  computational  labor  is 
enormous.   Modern  computing  equipment  is  helping  break  down  the  second  barrier 
and  we  hope  that  this  work  will  help  break  down  the  first. 

2.   PATH  INTEGRAL  FORMULATION  FOR  THE  SLATER  SUM 

The  path  integrals,  or  more  explicitly,  the  conditional  Wiener 
integrals  in  terms  of  which  we  express  the  Slater  sum,  may  be  defined  in 
the  following  way.   Let  the  parameter  t  be  defined  on  the  interval  (0,  p) 
and  let 

r(x)  =  (xx(t),  y^T),  zx(t),  ..,,  xn(t),  Yn(t),  zn(t))        (2.1) 

denote  a  continuous  function  of  T,  with  the  condition  r(0)  =  0:   it  is 
convenient  to  picture  rff)  as  the  generator  of  a  path  in  the  3N-dimensional 
coordinate  space  of  the  system  as  T  goes  from  0  to  p.   Let  F[r(x)]  denote  a 
functional  of  r(i).   Finally,  let 

r(x;n)  =  (x^Tjn),  y^i^n),    z^Tjn),  ...,  xN(x;n),  yN(-r;n),  zN(x;n))   (2.2) 
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denote  an  r(x)  which  is  piecewise  straight  and  has  breaks  at 

T  ,    1   ,    ...,    i        ;    such  a  function  is  displayed  in  Fig.  1  for  n  =  k   and  one 

space  coordinate.   The  end-points  of  the  T-interval  are 

T0  =  °'    Tn  =  P-  (2,3) 

The  conditional  Wiener  integral,  E{F|r(p)  =  R},  of  the  functional  F[r(T)] 
is  defined  by 

+oo       4-oo 

E{F|r(p)  =  R)  =  lim  J    ...   /  F[r(x;n)  ]d^,  (2.4) 


n-*» 

-00  -00 


where 


du  =A  n  J(2*(T.  -,  -T.))  2  exp    ~1+1   gV  r  n  d?-      (2-5) 

n    ni=0  ^    1+1   X  2(Ti+l  "  Ti}  V    1=1   ~x 


(r.  ,  -  r.)  =  Z^(x.(t.  ,;n)  -  x.(T.;n))  +  (y.(ir.  ,  :n)  -  y.(T.;n)) 
^i+l   ~i7    -_i  I  J  1+1      J  1         J  1+1      J  x 

J 

(2.6) 


+  (^/V'n)  -  zj<Vn))  J  ; 


dJ  r  =  n  dx  (t  ;n)  dy  (t  ;n)  dz  (t  ;n);  (2.7) 

--lj-1-      J-1-      J-L 


and  A  is  chosen  so  that 
n 


4-00        4-00 


du-  =  1.  (2.8) 

n 
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Piecewise  continuity  of  F[r(i)]  is  sufficient  to  ensure  the  existence  of 

the  limit  in  Eq.  (2.4),  and  we  will  assume  that  this  sufficiency  condition 

is  satisfied. 

We  wish  to  draw  attention  to  two  features  exhibited  in  the  above 

relations.   The  measure  du  is  a  Gaussian  probability  and  so  one  may  loosely 

regard  the  conditional  Wiener  integral  as  an  average  of  the  functional  F, 

where  the  average  is  taken  over  all  paths  r(f),  with  the  property  r(0)  =  0, 

r(p)  =  R:  the  statistical  weight  of  a  path  being  characterized  by  the  fact 

that  the  infinitesimal  increments  r(T  +  ot)  -  r(f)  are  governed  by  a  Gaussian 

distribution.   The  measure  du  is  invariant  to  the  transformation 

n 

T'  -   T,  r'(T')  -  r(T),  (2.9) 

where 

r(T)  =  \Zar'(T'  ),  (2.10) 

t  =  a  t1 .  (2.11) 

Hence  the  T-interval  can  always  be  normalized  to  (0,  l)  by  a  change  in  the 
space  coordinates . 

The  basic  relation  connecting  the  conditional  Wiener  integral  to 

3 
the  Slater  sum  is 


-f  (    (H«  -  *  A         f 

(2*0)  '   exp   -  ~  gp  ~   JE(exp  (- J   Vn(p(t)  +  R)dT)|r(p)  =  R'  -  R} 

V  /         0 

(2.12) 

=  Zf.  (R)Y.(R')  e   1, 
i 


where  on  the  left  the  T-interval  is  (0,(3),  R,  R'  are  two  points  in  the 


*•*»     *^t 


2 
3N-dimensional  coordinate  space,  and  (R*  -  R)   is  the  square  of  the  distance 

between  them  (cf.  Eq.  (2.6));  and  on  the  right  the  sum  extends  over  all 

eigenstates,  characterized  by  eigenvectors  ¥.  and  eigenvalues  E.  of  the 

equation 


|  (  S  ^2   +  ^2   +  ^2  )  *  +  (E  "  V  ¥  =  °'  (2'13) 

\i=l  dx.     dy.     oV 

where  V  is  required  to  give  a  discrete  spectrum. 

Equation  (2.13)  is  the  Schrodinger  equation  in  units  chosen  to 
make  m  =  H  =  1,  and  we  now  adopt  these  units;  the  thermal  wave  length  is  now 
v2np.   It  is  to  be  noted  that  the  potential  energy,  V  ,    appears  in  the 
functional 


F[r(T)]  =  exp  (-  /  V  (r(x)  +  R)cLt)),  (2.14) 


which  is  averaged  over  all  paths.   Because  r(i)  +  R  appears  in  the  argument 
of  V  ,   we  may  picture  the  effective  path  as  one  in  which  the  system  goes  from 
point  R  to  point  R'  in  the  T-interval  (0,(3).   This  picture  appears  again 
when  it  is  recognized  that  the  two  sides  of  Eq.  (2.12)  are  two  ways  of  writing 
the  Green's  function  for  the  Bloch  equation 


".♦-I-  (2-15> 
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The  above  formulation  ignores  symmetry  conditions  on  the 
eigenfunctions .   Taking  these  conditions  into  consideration,  let  P  denote 
the  permutation  operator,  and  construct  symmetric   and  antisymmetric 
eigenfunctions , 


1     .  =  -=■  £  P  ¥.(R),  (2.16) 


S'J  n/n:  P   J 


l.-pEopP  y.(g),  (2.17) 

where  a  =  +1  or  -1  according  as  the  permutation  is  even  or  odd.   Now  by 
applying  P  to  both  sides  of  Eq.  (2.12),  with  the  convention  that  P  operates 
on  the  primed  coordinates,  and  summing  over  P,  one  obtains 


1 


-3|       (RI  _  R)2  v         £ 
Z  P(2*0)     exp(-   ~    -   JE(exp  (-  /  VH(r(«r)  +  R)cLt)  |r(p)  =  R'  -  R) 


JWl     P  -      V     2P 


*  -PE. 

=  Z  y  (r)  y  e   J,  (2.18) 

J 


-f    /   (H-  -  R)2|       /? 
—  Z  opP(2rf)    exp(-   ~    ~  JE(exp(-  /  VN(r(-r)  +  g)dx)|r(p)  =  R"  -  R) 

"v  N .   P  _ 


*  -PE. 

=  ZY.  (r)Y   .  e    J.  (2.19) 

,   J  ~  a, j 
J 
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We  now  impose  the  requirement  that  R  be  a  point  derived  from  R  by  some 

permutation  of  the  particle  coordinates.   Noting  that  for  any  E.  there  are 

N!  Y.'s  because  of  the  symmetry  degeneracy,  and  that  in  a  sum  over  these 
J 

¥.'s,  ¥   .is  constant  and  ¥   .  changes  sign  we  have  the  fundamental 
relations 

5'  "  S)2  f  ,  b 

Z  P  exp  ( gg )E{exp(-  J   VN(r(T)  +  R)dT)|r(fs)  =  R  -  H)  =  V^,    (2.20) 

P  0 

(s'  -  s?2  § 


£  ap  P  exp(-  -  2p  -   )E(exp(- J   Vn(p(t)  +  g)dT)|r(0)  =  g'  -  R)  =  wj,   (2.21) 
P  0 


"R         F 
where  W7I  and  W   are  the  Slater  sums  for  Bose-Einstein  and  Fermi-Dirac 

statistics.   The  relation  shown  in  Eq.  (2.20)  was  used  some  years  ago  in  a 

study  of  the  ^-transition  of  helium. 


3-   FORMULAS  FOR  THE  CALCULATION  OF  W0 

2 

To  compute  W_  it  is  convenient  to  use  the  center  of  mass  and  relative 
coordinates.   The  Slater  sum  separates  into  the  product  of  two  Slater  sums,  one 
for  the  center  of  mass  coordinate  R,  which  can  be  evaluated  immediately,  and 
one  for  the  relative  coordinate  S, 


-11- 


k2 

s       n      j 3,         -ik  •    R        -p  f"     ik  •  R  „  -pE. 

W     =2A.6      /     ^e     ~       ~e  4e~  ~     Z  0*(s)   0. (s)    e        \ 

2  J       (2*)3  i      X  ~ 


22  \3Z0.    (S)    0.(s)    e        \  (3-D 

i 


where  0.(s)  satisfies  the  Schrodinger  equation  for  the  relative  coordinate, 


V20.(S)  +  (E   -  V(S))  0.(S)  =  0.  (3-2) 

l  ~      l     ~    l  ~ 


The  interchange  of  the  two  particles  does  not  alter  the  center  of  mass 
coordinate  but  changes  the  relative  coordinate  S  into  -  S,  so 


Wg=W°  +  ^,  (3-3) 


w^  =  w°  -  w^,  (3.4) 


where  the  direct  term,  W  ,    is  given  by 


3 

W°  =  22  K3     E  0*(s){Zf.  (S)  e   \  (3-5) 

£.  l  ~  l  ~ 

l 


and  the  exchange  term  WT  is 


3 

—  —  P)V 

W*  =  22*.3Z0*(-S)  0.(s)  e    \  (3-6) 
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It  is  to  "be  recognized  that  the  sums  in  Eqs .  (3-5)  and  (3-6)  extend  over 
all  states  without  regard  to  symmetry.   Incorporating  the  physical  constants 
into  Eq.  (2.12)  and  performing  the  normalization  of  t  to  the  interval  (0,1) 
according  to  Eqs.  (2-9),  (2.10),  and  (2.1l)  one  obtains, 


W°  =  1 


|  exp  f-  p/v(p  r(T)  +  S)dTJ  |r(l)  =  Cij,  (3-7) 


and 


2ttS2 

TTE       \2 
W2  =  e 


1 
EJexp  (-P/V(j"£(^)  +  S)dTJ|r(l)  =  -^  s}  ,     (3-8) 


where  r(x)  is  the  generator  of  the  path  of  the  system  in  the  3 ""dimensional 
relative  coordinate  space  as  T  goes  from  0  to  1.   A  transformation  can  be 
performed  on  r(x)  in  Eq.  (3«8)  "to  make  the  condition  at  r(l)  the  same  as 
that  in  Eq.  (3-7)-   The  result  is 


2nS2 

,2 
E       ^ 
W*  =  e      E 


J  exp  (-pj  V(jz   t(t)  +  S  -  2TS)dTJ|r(l)  =  o|  .    (3-9) 


Let  F[r(f)]  denote  either  of  the  two  functionals, 


1  \ 


F°[r(T)]  =  exp   -p  /  V(^r  r(x)  +  S)dT*  ,  (3-10) 


0 


or 


■13- 


FE[r(T)]   =  exp  (-p/v(^r(T)  +  S 


2TS)dT   .         (3-11) 


Then  the  numerical  scheme  for  evaluating  the  conditional  Wiener  integrals 
in  Eqs  .  (3*7)  and  (3-9)  is  to  approximate  each  "by  a  3n-dimensional  integral, 


E[F|£(1)  =  0}  ZJ    ...  J   F[r(T;n)]dun,  (3-12) 


where  du,  is  given  in  Eq.  (2.5)  with  r_  =  r  =0.   The  break  points  in  the 
n  ~0   ~n 

piecewise  straight  path  r(f;n)  are  chosen  at  equal  time  intervals,  T.  =  — , 

i  =  0,  1,  .  ..,  n.   The  3n_dimensional  integral  is  then  evaluated  by  a  Monte 

Carlo  sampling  procedure.   The  sampling  is  done  by  choosing  the  coordinates 

of  the  break  points  r.  =  r(f.;n)  of  r(x;n)  according  to  the  distribution 

dp.  .   F[r(i;n)]  is  then  evaluated  with  this  path  for  a  set  of  values  of 

S  =  |S|.    Further  piecewise  straight  paths  are  then  chosen  and  F[r(T;n)] 

is  averaged  over  all  paths  for  each  value  of  S. 

Concerning  the  choice  of  the  r.,  it  is  evident  that  the 

~i 

distribution  du.  does  not  give  independent  Gaussian  increments  (r.  ,  -  r.  ) 
n  B  v~i+l   ~iy 

due  to  the  condition  r  =  0:  however,  the  choice  of  the  r.  can  be  made  to 

~n  ~i 

depend  on  independent  Gaussian  random  variables  by  use  of  an  interpolation 

6 
formula  for  a  conditional  Brownian  motion  path.     If  r_,  r,  ,  •  ..,  r.  are 

~0  ~1       ~i 

fixed  then, 


r.  (t   -  t.  ,)  +  r   (t.  ,  -  t.)  /(t.  .  -  t.)(t  -  T.  ,) 

~i  v  n    l+l '       ~n  v  l+l    \J  _/v  i+I    \JK   n    i+l' 

r.  ,  =  +  ?  Y 

~1+1                T   -  T .  <-   '          T   -  T . 


n    i  n    i 

(3-13) 


■Ik- 


where  the  coordinate  random  variables  of  £  =  (£  ,    £  ,  £  )  are  independent 
and  Gaussianly  distributed  with  mean  0  and  variance  1. 

The  labor  involved  in  the  evaluation  of  F[r(i;n)]  is  reduced  by 

the  fact  that  it  is  not  necessary  to  evaluate  the  T-integral  to  a  high 

7 
order  of  accuracy.   On  the  basis  of  some  earlier  work  we  can  expect  the 

approximation  represented  by  Eq.  (3-12)  to  have  an  error  not  less  than 


m 


8 


It  is  therefore  sufficient  to  evaluate  the  T-integral  by  applying 


the  trapezoidal  rule  to  the  intervals  (t.    t.),  i  =  1,  2,    ...,    n.   Then, 


FD[r(T;n)]sexp  [-£  "z"  V(£  r. 

i=0   >/« 
\ 


+  s)  • 


(3.1*0 


Let  M  paths  r  (i;n),  j   =  1,  2,    . . . ,   M  be  chosen  according  to  the 
sampling  scheme  described  above,  let  V  be  given  by  Eq.  (l.5)>  and  introduce 
the  dimensionless  variables, 


/Jit 


(3-15) 


D       F 
Then  the  numerical  approximations  to  W  and  W  are  given  explicitly  by  the 

formulas , 


M 
W2"  M  h     exp 


n-1 

*  z 

n  i=0 


T3.     +    ( 

1 


•12 


A 


(3.16) 


and 
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M 
W2  ~  M     L     6XP 


r1?   +  d(l-2T.  ) 


■12 


'rJ   +  d(l-2T.  ) 


(3-17) 


4.      RESULTS   OF  THE   CALCULATION  OF  W 


B 


It  is  evident  from  the  formulas  in  the  last  section  that  they 
approach  the  classical  result  when  the  thermal  wave  length  vanishes;  i.e., 
when  \   ->  0,  Eqs .  (3-7)  and  (3*8)  give 


W 


D 


exp  (-|3V(S)), 


(k.l) 


wJ-0.  (4.2) 


One  can  see  from  a  simple  qualitative  argument  that  for  fixed  n  the 
approximation,  Eq.  (3-12),  is  expected  to  improve  as  \   ->  0.   In  this 
approximation  the  true  ensemble  of  paths  has  been  replaced  by  an  ensemble 
of  broken  straight  line  paths,  one  straight  line  segment  in  such  a  path 
corresponding  to  a  "time"  interval  of  —  in  length.   If  we  consider  a  finer 
subdivision,  obtained  by  chopping  each  of  the  n  intervals  in  half,  then 
each  straight  line  segment  will  have  a  break  at  the  center  as  illustrated 
in  Fig.  2  for  one  space  coordinate.   We  can  now  think  of  the  deviation,  |, 
at  the  center  of  the  first  interval  in  Fig.  2  as  a  random  variable.   From 
the  interpolation  formula,  Eq.  (3-13);  it  follows  that  it  may  be  regarded 
as  a  Gaussian  random  variable  with  variance  i — .   Thus  in  our  original 
approximation,  with  time  segments  — ,  we  are,  roughly  speaking,  ignoring 
fluctuations  in  each  space  coordinate  of  the  order 
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K_ 

v  it 

\ 

the  factor  enters  because  it  multiplies  the  path  coordinate  in 

v  n 
Eqs  .  (3-7)  and  (3-8).   It  is  therefore  reasonable  to  assume  that  the 

accuracy  of  our  approximation  will  improve  as  \   -»  0.   Conversely,  we 

can  expect  a  larger  error  as  \   -»  °°;  i.e.,  as  the  temperature  becomes 

small.   The  above  argument  suggests  that  a  decrease  in  temperature  must 

be  compensated  by  an  increase  in  n  such  that  Tn  remains  constant,  if  the 

error  is  to  remain  constant  as  the  temperature  is  lowered.   Since  the 

computing  time  depends  critically  on  n,  it  is  clear  that  one  cannot 

expect  to  pursue  these  calculations  to  arbitrarily  low  temperatures. 

Computing  time  restrictions  led  us  to  use  the  value 

n  =  512  =  29  (k.3) 

in  almost  all  of  the  calculations.   It  is  extremely  difficult  to  get  an 
a  priori  estimate  of  the  error  to  be  expected  for  a  given  temperature  and 
value  of  n.   However,  one  can  get  a  useful  picture  of  the  error  in  the 
following  way.   As  a  consequence  of  the  argument  in  the  last  paragraph  we 
certainly  must  restrict  our  attention  to  temperatures  such  that  fluctuations 
in  position  of  the  order 


\        -y  71 


n 


2~5^3  (k.k) 


-17- 


are  small  compared,  with  the  range  of  the  potential.   Retaining  0.1$ 
(relative  to  the  potential  minimum)  accuracy  in  the  potential  then,  from 
Eq.  (1.9),  it  is  reasonable  to  regard  the  range  of  V  as  ha  ;    i.e.,  to 
regard  V  -   0   for  r  >  ha.      Hence  our  criterion  becomes: 


0 


must  be  small,  relative  to  unity,  or,  substituting  the  values  of  k  and  o, 

— '- must  be  small  relative  to  unity.   This  crude  argument  ignores  the 

n/~T 

r 
fact  that  V  changes  very  rapidly  when  —  <  1,  however  this  neglect  is  not 

-BV 

so  serious  as  it  might  appear  since  e    is  practically  zero  for  these 

values  of  r. 

E 

Except  at  low  temperatures,  the  exchange  term  W~  should  be  small. 

An  estimate  of  the  upper  bound  for  this  term  is  easy  to  construct  from  the 
present  formulation.   To  get  this  estimate  we  make  a  slight  change  in  the 
potential  and  insist  that  V(r)  =  °°  for  r  <  0  and  otherwise  is  given  by 
Eq.  (1.5).   Under  this  condition  the  Wiener  integral  factor  in  Eq.  (3-9) 

is  zero  for  |S|  <  a   and  cannot  exceed  e   for  |S|  >  a.   It  follows  that 

E  Q 

Wp,  for  this  modified  potential,  satisfies  the  inequality 

2™    +  pa 


W?<e    X  ,  (i+.5) 


2 


and  substitution  of  the  values  of  the  physical  constants  gives 
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-0-5^T  +  10 
W*<e  T  (4.6) 

E 

Our  results  show  that  W  is  considerably  smaller  than  this  bound.   The 

most  likely  explanation  for  this  is  the  following  one.   Examination  of 
Eq-  (3-9)  shows  that  at  the  midpoint  of  the  path  (t  =  — )  the  potential  is 


T<7=-s(i». 

v  it 


therefore  it  can  be  expected  that  the  Wiener  integral  in  Eq.  (3-9)  will  be 
very  small  except  when  \   is  large  enough  to  make 


7-  r2(|)  >  a2  (4.7) 

2  1 
with  a  nonvanishing  probability.   If  we  simply  replace  r  (— )  by  its  mean 

o 

value,  j-   (see  Eq.  (3*13) )>  this  inequality  becomes 

T  <  2.8°  K.  (4.8) 


This  argument  suggests  that  we  can  only  expect  a  significant  contribution 
from  the  exchange  term  below  2.8  K,    and  our  numerical  results  support 
this  conclusion.   This  argument  points  to  an  interesting  picture  of 
contributions  to  the  exchange  term.   One  contribution  comes  from  the 
exponential  factor  which,  as  we  have  seen,  introduces  a  factor  of 
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2ng2 
e 


in  the  bound.   Thus  this  factor  will  make  the  exchange  term  go  to  zero  like 

-kT 
e    as  T  -*  0°.   One  can  picture  this  contribution  as  coming  from  the  kinetic 

energy  associated  with  the  relative  motion  when  the  particles  exchange 

position.   Another  contribution  comes  from  the  interaction  during  the  exchange 

in  position  and  this  contribution  also  makes  the  exchange  term  vanish  as 

T  -»  oo ;  see  Fig.  3-   Our  results  suggest  that  the  latter  effect  is  probably 

more  important  than  the  kinetic  energy  effect  in  making  the  exchange  term 

small . 

In  Fig.  h   the  Slater  sum,  W  ,  as  a  function  of  —  for  different 
temperatures  is  displayed.   Some  results  have  been  omitted  to  improve  the 
legibility  of  this  figure.   Only  the  curve  for  T  =  2  K  explicitly  includes 
the  exchange  term.   At  T  =  5  K  the  exchange  term  was  found  to  be  negligible 
compared  with  the  direct  term  so  it  was  not  calculated  for  the  higher 
temperatures  and  is  omitted  in  the  results  for  T  =  5  K  and  above.   The 
location  and  height  of  the  maximum  is  given  in  Table  I. 

There  are  two  important  sources  of  error  in  these  calculations ,    one 
arising  from  the  Monte  Carlo  sampling  and  the  other  arising  from  the  approxi- 
mation represented  in  Eq.  (3-12).   For  reasons  just  discussed  it  is  to  be 
expected  that  these  errors  will  be  most  important  at  lower  temperatures.   A 
measure  of  the  Monte  Carlo  sampling  error  is  shown  in  Figs.  5a  and  b  where 
the  calculated  points  with  standard  deviations  are  shown  for  T  =  5  K  and 
30  K.   The  length  of  the  vertical  line  at  a  point  is  twice  the  standard 
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deviation  of  the  sample  consisting  of  1000  independent  paths.   Standard 
deviations  are  not  shown  in  the  tails  where  they  are  too  small  to  draw 
on  this  scale.   Some  measure  of  the  other  error  is  given  by  comparing 
results  based  on  other  values  of  n,  the  number  of  straight  line  segments 
in  the  path.   In  Fig.  6  the  results  of  calculating  Wp  (the  direct  term) 

at  T  =  2°  K  for  n  =  100  and  n  =  512  are  shown. 

^  o 

In  our  initial  calculations  we  only  went  down  to  T  =  5  K 

because  we  suspected  that  the  errors  encountered  at  lower  temperatures 

would  be  too  large.   After  learning   that  Larsen  and  Kilpatrick  had 

calculated  W.  at  T  =  2  K  in  an  entirely  different  way  we  were  stimulated 

to  push  our  calculations  down  to  2  K.   The  results  compared  with  those  of 

Larsen  and  Kilpatrick  are  shown  in  Figs.  7a  and  b;  the  direct  term  is 

shown  in  Fig.  7a  and  the  exchange  term  is  shown  in  Fig.  7b.   The  sampling 

error  at  a  representative  set  of  points  is  shown. 

It  is  of  some  interest  to  compare  these  results  against  the  pure 

classical  value  of  W  and  its  value  from  the  first  few  terms  of  an  expansion 

in  powers  of  X,    sometimes  called  the  Wigner-Kirkwood  expansion.   The 

classical  value  is  given  by 


W2(classical)  =  e"^V.  (^.9) 


The  Wigner-Kirkwood  value  displayed  here  is  obtained  by  truncating  the 
expansion  at  terms  in  —-y-     to  obtain 


W  (W-K)  =  W  (classical)   •  C,  (^-10) 


11 
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where 

C  =  exp  (2.5  \   (f)8  ♦  7.0  \   (f)10  ♦  (1.5  %  "  U  %)(f)^)^      C*-H) 

P  P  P  P 

and  where  7  and  p  are  given  in  Eq.  (3-15) •   This  comparison  is  made  in 
Figs.  8a,b,c  at  T  =  2  ,  10  ,  and  30  K.   Only  the  direct  term  is  shown. 

One  cannot  of  course  expect  good  agreement  but  it  is  nevertheless 
tempting  to  compare  the  radial  distribution  function  in  the  present 
approximation  against  experimental  results.   This  comparison  is  made  in 
Table  II  where  we  compare  location  and  height  of  the  first  maximum  in  g(s) 

and  the  point  at  which  g(s)  first  increases  above  zero  for  He,  according 

12 
to  measurements  by  Henshaw   with  our  results.   We  note  that  although  there 

is  a  marked  difference  in  the  height  of  the  maximum,  the  other  parameters 

are  in  fairly  good  agreement.   It  is  to  be  noted  that  the  position  of  the 

maximum  moves  to  higher  S  as  the  temperature  increases  while  our  results 

show  the  opposite  behavior.   This  qualitative  difference  is  almost  certainly 

due  to  the  fact  that  the  approximation  used  here,  Eq.  (1.19),  includes  no 

dependence  on  density.   Henshaw' s  measurements  were  taken  at  T  =  2.2  K, 

density  =  0.1U6  gm/cirr   and  T  =   5-0^°  K,    density  =  0.095  gm/cm    . 

The  second  virial  coefficient,  Eq.  (l.20),  obtained  from  our 

calculations  is  compared  with  results  by  Kilpatrick,  Keller,  Hammel  and 

13 
Metropolis    in  Table  III.   The  results  marked  by  an  e  were  obtained 

Ik 

by  us  using  the  high  temperature  expansion.    Experimental  results  are 

displayed  in  column  k. 
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T(°  K)  Max  W? 


S 

2  o 

2                   1.94  1.52 

5                 1.55  1.^6 

10                 1 .  ko  1 .  39 

20                 1.28  1.31 

30                    1.22  1.27 

kO                                                  1.18  1.24 

50                    1.15  1.22 

75                    1.11  1.19 

100                                      1 . 09  1 . 18 

273.18                              1.04  1.11+ 

"P 

Table  I:   Maximum  value  of  W   (column  2)  and  the 

location  of  this  maximum  (column  3)  for 

different  temperatures  (column  l). 
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T(° 

K) 

This  Work 

2 

-182.84 

5 

-  59-41 

10 

-  21.30 

20 

-  2.49 

30 

3.68 

40 

6.59 

50 

8.26 

75 

10.28 

100 

11.08 

273- 

18 

11.65 

B 


Kilpatrick  et  al 

-177.39 

-  59-14 

-  21.34 

-  2.53 
3-57 
6.49 
8.16 

10 .  l4e 
11 . 02e 
11.596 


Experiment 


-193.3 


62.2 
23.4 


4.04 


b 


2.42  D 

6.57(4o.09°k)c 

8.o6(50.09°k)c 

10. 70(75. oi°k)c 

11.85(100. 02°K)C 
11.77  d 


Table  III:   The  second  virial  coefficient  compared  with  the  computer  results 
of  Kilpatrick,  Keller,  Hammel  and  Metropolis  and  with  experiment. 

a.  Obtained  by  linear  extrapolation  of  the  data  in  W.  E.  Keller, 
Phys.  Rev.  97,  1  (1955). 

b.  David  White,  Thor  Rubin,  Paul  Camky  and  H.  L.  Johnson, 
J.  Phys.  Chem.,  64,  1607  (i960). 

c.  W.  H.  Keesom,  Helium   (Elsevier,  Amsterdam  1942). 

d.  W.  G.  Schneider  and  J.  A.  H.  Duffie,  J.  Chem.  Phys.  17, 

751  (1949). 

e.  Obtained  by  us  from  the  high  temperature  expansion 
(footnote  reference  l4). 
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FIGURE  CAPTIONS 

Figure  1:    Example  of  x(i;k). 

Figure  2:    Illustration  of  the  error  in  using  paths  x(t;^).   Fluctuations 
represented,  by  £  are  ignored. 

Figure   3:    Two-dimensional  illustration  of  the  qualitative  difference 

between  a  high  temperature  exchange  path  and  a  low  temperature 
exchange  path  (relative  coordinates).   As  the  temperature  becomes 
higher  the  path  tends  to  follow  more  closely  the  straight  line 
connecting  points  S  and  -S. 

Figure  h:      W^asa  function  of  -  for  T  =  2°  K,  5°  K,  10°  K,  75°  K 
and  273.18°  K. 

Figure  5a:   W  as  a  function  of  —  for  T  =  5  K  with  the  Monte  Carlo  sampling 
error  shown.   Total  length  of  the  error  indicator,  represented  by 
the  vertical  line  segment,  is  twice  the  standard  deviation  of  the 

sample  mean,  represented  by  the  dot. 

"R  *-?  n 

Figure  5b:   W  as  a  function  of  —for  T  =  30  K  with  the  Monte  Carlo  sampling 

error  shown.   Total  length  of  the  error  indicator,  represented 

by  the  vertical  line  segment,   is  twice  the  standard  deviation  of 

the  sample  mean,  represented  by  the  dot. 
Figure  6:   W   (the  direct  term)  at  T  =  2°  K  for  n  =  100  and  n  =  512. 
Figure  J&:        Comparison  of  our  results  with  those  of  Larsen  and  Kilpatrick 

for  the  direct  term  at  T  =  2  K.   The  Monte  Carlo  sampling  error 

(as  in  Figs.  5a, b)  is  shown. 
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FIGURE  CAPTIONS  (CONT'D) 


Figure  7b:    Comparison  of  our  results  with  those  of  Larsen  and  Kilpatrick 
for  the  exchange  term  at  T  =  2  K.   The  Monte  Carlo  sampling 
error  (as  in  Figs.  5a, b)  at  a  representative  set  of  points  is 
shown. 

Figure  8a:    Comparison  of  the  classical  approximation,  a  Wigner-Kirkwood 
approximation  (Eq.  4.10),  and  our  results  for  W  at  T  =  2  K. 

Figure  8b:    Comparison  of  the  classical  approximation,  a  Wigner-Kirkwood 

approximation  (Eq.  4.10),  and  our  results  for  W  at  T  =  10  K. 

Figure  8c:    Comparison  of  the  classical  approximation,  a  Wigner-Kirkwood 
approximation  (Eq.  4.10),  and  our  results  for  W  at  T  =  30  K. 
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Figure  2. 
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